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ABSTRACT 

The stability theory for finite difference Initial Boundary- Value approximations to sys- 
tems of hyperbolic partial differential equations states that the exclusion of eigenvalues and 
generalized eigenvalues is a sufficient condition for stability. The theory, however, does not 
discuss the nature of numerical approximations in the presence of such eigenvalues. 

In fact, as was shown previously [1], for the problem of vortex shedding by a 2-D cylinder 
in subsonic flow, stating boundary conditions in terms of the primitive (non-characteristic) 
variables may lead to such eigenvalues, causing perturbations that decay slowly in space and 
remain periodic time. Characteristic formulation of the boundary conditions avoided this 
problem. 

In this paper, we report on a more systematic study of the behavior of the (linearized) 
one-dimensional gas dynamic equations under various sets of oscillation-inducing “legal” 
boundary conditions. 
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1. Introduction 


The increase in computers’ speed and memory allows researchers to investigate fluid 
dynamical problems with greater attention to delicate features of the flows. For example, 
only in the last five years have investigators [1,6] computed the vortex shedding phenomenon 
behind a two-dimensional cylinder at low Mach numbers and moderate Reynolds numbers. 
Unfortunately, with this improved ability due to increased computer power to compute 
complex phenomena comes the necessity to deal with unforeseen numerical surprises which 
might be mistaken for real or physical effects. For example, in the case of the 2-D cylindrical 
Von-Karman vortex street computation - the main features of the flow, such as the shedding 
frequency are predicted accurately. However, a concomitant computational result appearing 
in the output data base is a spurious secondary frequency which was mistakingly attributed 
to the start of transition. In previous work, [1], we have shown that the spurious secondary 
frequency resulted from applying the far-field inflow boundary condition to the primitive 
variables. It was also shown there that the same boundary treatment, applied, however, 
to the characteristic variables, eliminated this phenomenon. In this paper, we characterize, 
by an analytic description, the boundary conditions (both inflow and outflow) under which 
the numerical solution of the Euler equations will exhibit temporal oscillations which are 
foreign to the exact solution of the p.d.e.’s. The starting point of our study is the modal 
analysis developed by G-K-S [2], and Osher [3]. This theory states that for finite differences 
approximations to Initial Boundary Value problems of hyperbolic systems stability is assured 
by the exclusion of eigenvalues and generalized eigenvalues. The theory, however, does not 
discuss the nature of the numerical approximations in the presence of such eigenvalues. 

The model that we study is that of the 1-D compressible Euler equations linearized about 
free stream conditions. The numerical schemes that we analyze are second order in space 
and time represented by the Lax-Wendroff scheme. In our 1-D case it is equivalent to all 
other second order algorithms such as the MacCormack scheme [4]. 

In Section 2, we briefly review the underlying theoretical considerations both for the 
p.d.e. and the f.d.e. formulation. 

In Section 3, we discuss the case of inflow b.c.’s expressed with primitive variables. We 
show that even though the numerical solution is technically stable and therefore, by the Lax 
equivalence theorem, convergent, for finite meshes and time it shows spurious oscillations 
that decay only slowly with mesh refinement. In Section 4, we repeat this demonstration 
for the case of primitive- variable description of the out-flow boundary conditions. In Section 
5, we apply a primitive- variable formulation of the boundary conditions both at the inflow 
and outflow boundaries. This allows a nonlinear interaction between the modes created at 
both boundaries. At this point, the reader should be reminded that the same treatment 
of the boundary conditions, applied to the characteristic variables eliminates any spurious 
frequencies. This is a corollary of our analytical formulation and is borne out in our numerical 
experiments. 

All sections contain numerical examples pertinent to the initial boundary value problem 
considered therein. 
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2. Analytical Preliminaries 


We start with the one- dimensional Euler equations of gas dynamics in conservation form: 


where 


dU dF(U) 
dt dx 


= 0 


U = 


P 

m 

E 


( 2 . 1 ) 


( 2 . 2 ) 


F = 


m 



u(E + p) 


(2.3) 


Here u,p,p and E are, respectively, the velocity, density, pressure, and total energy per unit 
volume, m = pu is the mass flux. In the case of an ideal gas the equation of state is 


P- (7-1) 




(2.4) 


where 7 = 7/5 for diatomic gases such as air. Equation (2.1) can be rewritten in non- 
conservation form as 

§ + ^)f = 0 (2.5) 

where A = dF/dU is the Jacobian of the flux vector F with respect to the solution vector 
U. Linearizing about steady free stream conditions, = (/>«>, PooV-oo, -Eoo)> Equation (2.5) 
becomes: 

gj(«0 + A(U^(SU) = 0, (2.6) 

where 6 U = U — I/qo is the perturbation vector. The matrix A(Uoo) = [ ^f^ ] has 
three eigenvalues, — c^, a 2 = -f and a 3 = The free stream sound 

speed, Coo, is given by Coo = (7P00/ Poo) 1 • The corresponding eigenvectors in terms of the 
conserved-variable perturbations are given by 


Ri = — ( 8m — u^p ) + 


7-1 


^ u lo S P - «oo Sm + 6e] 


R 2 = (6m — UooSp ) + 


Coo L2 

7 1 ^ u lof>P ~ Uoo6m + 6E 


and 


R 3 = Cqq 6p - 


7-1 


L2 


-u^p — Uoo6m + 6E 
• £ 


Furthermore, using a linearized version of the equation of state, i.e,, 

n 


= (7 - 1) 


- u^Sm + SE 


(2.7) 

( 2 . 8 ) 

(2.9) 


( 2 . 10 ) 
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equations (2.7) - (2.9) may also be written in terms of the primitive- variable perturbations, 
(6p,8u,6p): 

R\ = —PooSu+ — 8p , ( 2 . 11 ) 


and 


R 2 = Poo8u -) 8p, 

Cqo 


R 3 = Coo 8p — —8p. 


(2.12) 


(2.13) 


Any perturbation imposed on the free stream solution will evolve as a combination of these 
eigenvectors. 

In terms of the characteristic variables R,(s = 1, 2, 3) defined in (2.7 - 2.9) and/or (2.11 
- 2.13), Equation (2.6) may be written as: 


or 


dRi . .dR x 

-ar + (u ~- c ~>-& =0 

dR 2 , 

dt + ^°° + Coo) dx ~ ° 

SR* . 9R 3 

-er + ““-fl7 = ° 

SR, SR. ■ 

+ Q.-zt = 0, a = 1,2,3. 


(2.14) 

(2.15) 

(2.16) 

at ' '• ax ^ 2 ' 17 ' 

In a finite domain, say 0 < x < 1, for the subsonic case u M < c^,, the system (2.7) is well 
posed with the following initial and boundary conditions: 


R a (x,0) = f,(x) s = 1,2,3 

R 2 ( 0, t ) - a o #i(0, t) = gi(t ) 

^3(0, t) — /3 o 72i(0,t) = g 2 (t) 

72i(l, i ) + ® 1 R % ( 1 1 i) + £i 7Z3(1) i) = 53(0 


(2.18) 

(2.19) 

( 2 . 20 ) 
( 2 . 21 ) 


where a 0 ,j3 0 , <j\ and t\ are arbitrary real constants and f t (x),g a (t ) (s = 1,2,3) are square 
integrable in their respective domains. 

We get numerical approximations of second order spatial and temporal accuracy by using 
the Lax-Wendroff scheme, 

< + ‘ = <-^(^-^,)+|^K«-2<+<,). i > 0 s= 1,2,3 (2.22) 


2(Ax) 


where w‘ ,n = w‘(jAx,nAt) is the finite difference approximation to R,(x,t). 

We note at this junction that although we shall illustrate the detailed development of the 
spurious frequencies using the Lax-Wendroff scheme, they depend in fact mostly on the form 
of the finite- difference boundary conditions (yet to be specified) and not on the particular 
inner algorithm. 
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To close the system (2.22) we need both the finite difference form of the algebraic bound- 


ary conditions (2.19 - 2.22) i.e., 

vjq™ — ct 0 wl’ n = <7i(nAt) (2.23) 

^o’" - /3 0 u4’ n = 92(,nAt) (2.24) 

w]f + oitojy" + £i w% n = g 3 (nAt ) (2.25) 

and also three “numerical” boundary conditions, one at x = 0 and two at x = 1 = AxN. A 
commonly used method of imposing “numerical” boundary condition is to extrapolate from 
the interior in the following manner: 

w 0 ' + co^o + £o^o = w i + &oV>i + £ 0^1 , (2.26) 

and at the outflow boundary, x = 1 = N Ax, 

w 2 ri n - aiw]f n = - attv/jp.!, (2.27) 

wff - = w]i n _ j - fawnU- (2.28) 


The “numerical” boundary conditions (2.26) - (2.28) are zeroth order extrapolations, with 
Ox , /?i , cr 0 and Sq arbitrary real constants. 

For stability analysis studies it is sufficient to consider the case of homogeneous boundary 
conditions, i.e., g r (nAt) = 0,r = 1,2,3. We look for solutions of the form 

w *•" = z n (A,K’ + B.ri), s = 1,2,3. (2.29) 

Substituting this ansatz into (2.22) we find that the k 3 ’s and /x 4 ’s are the roots of the 
quadratic equation 

(A? - \.)x] + 2(1 -z- \])x, + (A? + A.) = 0 (2.30) 

where X, = a,At/ Ax. For \z\ > 1 one of the roots of (2.30) (say « 4 ) is inside the unit circle 
and the other, fi 3 , is outside the unit circle. 

Substituting (2.29) into the homogeneous version of (2.23) - (2.28) we get the following 
system of equations: 

/ A x \ 

A 2 

Q(k>(z),M z )) £* =° (2.31) 

b\ 

\ b 3 ) 
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where Q(k,(z), n,(z))^Q(z) is given by 

«1 - 1 Co(k 2 ~ - 1 ) e o(«3 ” 1 ) /^1 — 1 <7 0 (/X 2 - 1) £ 0 (/i3 - 1) 

-Qo 1 0 — a 0 1 0 

-A) oi -p 0 oi 

GlK-2 fj-i crifJ-2 €1^3 

— C*1 «i r-1 («i - 1) «f _1 (K2 - 1) o - 1) - 1) 0 

-/3lKf _1 («l - 1) 0 K^ _1 (k 3 - 1) -Pitf-'im - 1) 0 - 1) 

(2.32) 

Note that the dependency on z in (2.32) comes from the fact that (j = 1,2,3) are 

functions of z defined implicitly in (2.30). 

It is quite easy to formulate a necessary condition for stability. This is the Riabenkii- 
Godunov condition [5]: 

Lemma 2.1: The Lax-Wendroff scheme (2.22) with the inflow boundary conditions (2.23), 
(2.24), and (2.26) and outflow conditions (2.25), (2.27), (2.28) is unstable if there is \z 0 \ > 1 
such that 

det Q(z 0 ) = 0. (2.33) 

In fact, such a z 0 gives a solution of the type 

w]' n = z^A,k\ + B.tf) 

that grows with the number of time steps . This is indeed a case of classical instability. 

The sufficient condition for stability can also be formulated in terms of the determinant 
of Q(z ); in fact it has been proven: 

Lemma (2.2): [GKS] The LW scheme (2.22) with the b.c.’s (2.23) - (2.28) is stable if for 
every |z 0 | > 1, 

det Q(z 0 ) -f- 0 (2.34) 

It is evident that the case |z 0 | = 1 with det Q(zq) = 0 is covered by neither lemma. It is 
precisely this case, however, that is responsible for the spurious oscillations. 

For the sake of convenience, we rewrite (2.31) and (2.32) in the following way: 
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Define 


/ «! - 1 a 0 (K 2 - 1) £o(*3 - 1) \ 


«w 


Now (2.31) becomes 


and 


A(k) = 


-a 0 

V -Po 


1 

0 


0 

1 


/ 


K 


N 


rr 

°1 K 2 


-a x «f 1 («1 - 1) K 2 X («2 - 1) 




V *(«! - 1 ) 



*(*3 - 1) y 



;4(re) o +A(/i) 6=0 
B(i c) a +i?(^) 6=0 


ow 


-A(«) ^ 

B(«) B^)) 


(2.35) 


(2.36o) 

(2.366) 

(2.37) 


The determinant condition (2.33) involves the solution of a very complicated nonlinear com- 
plex equation. Only little insight can be gained by trying to directly analyze it. One simple 
case, though, is given by the fully characteristic case 


CtQ — (3 q = (Jq — £q = Q-i — fix — Oi — £l — 0. 


(2.38) 


Note that in this case A and B are diagonal 



f «i - 1 0 0 ^ 


/ 0 0 \ 

A(k) = 

0 1 0 
1 0 oij 

£(*) = 

0 /C^ -1 («2 — 1) 0 

< 0 0 «f -1 («3 — 1) / 


and the determinant condition reduces to 

[(«i-l)/*f— (Mi — 1 )][(/*2 — l)/^ -1 — ^ _l ( K 2— 1)][(M3— — !)] 


(2.39) 


0. (2.40) 


Equation (2.40) indicates that the choice of parameters in (2.38) decouples the system (2.31) 
and the system is thus reduced to the scalar case. It may be easily shown, using (2.30), that 
there is no zq for which (2.40) has a solution, and therefore no spurious solutions exist in 
the fully characteristic case. 
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3. The Inflow Case 


In this section, we will investigate the possibility of having a spurious frequency induced 
by “primitive” boundary conditions (numerical and analytical) at the inflow boundary. By 
“primitive” boundary conditions, we mean b.c.’s expressed in terms of primitive variables. 
By a spurious frequency, we mean a solution of the type z 0 = e 1 * to the determinant condition 
(2.33). Indeed, if the determinant has such a root z 0 , then a spurious solution of the form 

< n = z n (A.Kl + B'A) = e tn *(A .> cj + B,fi) (3.1) 

is superimposed on any numerical solution of (2.22) - (2.28). 

When the number of grid point increases without limit, N — > oo, the influence of the 
outflow boundary treatment gets decoupled from that of the inflow boundary. This can be 
seen by observing that the submatrix B(i c) defined in (2.35) vanishes because |/cj < 1. At 
the same time B(fJ.) tends to infinity. Therefore, from (2.36b) we have that > 0. We have 
then, from (2.36a) 

A(k) a= 0 (3.2) 

and the determinant condition (2.32) reduces by (2.37) to 

det Q(z ) = (/ci - 1) + a 0 a 0 (K 2 - 1) + £oA )(«3 - 1) = 0. (3.3) 

In practice, one cannot achieve N —* oo. We then inquire how to best decouple the 
outflow influence from the inflow for the case of fixed N, albeit large. If we set the amplitude 
of the perturbations due to outflow treatment, | b | r equal to zero then from (2.36) we have 

(A(/c) + B(k)) a = 0 (2.36c) 

as well as A(k) a= 0. Thus the term B(i c) in (2.36c) represents the effect of the incomplete 
decoupling due to the finite N. Since N 1, it is clear from the definitions that B(k) 
represents a small perturbation to A(k), which we would like to minimize. It can be shown 
that this minimization is achieved for d = ai = ax = = 0. Notice that this corresponds 

to imposing both analytic and numerical outflow boundary conditions on the characteristic 
variables, R a (s = 1,2,3). 

A spurious solution exists if the determinant condition (3.3) has a root of the form 
z 0 = e“*\ In this case, the homogeneous version of (2.21) - (2.28) with oj = e 1 = = /3 X = 0 

admits a solution of the form 

w‘* = Ae‘”V„ ^0 

|K.(e~*)| < 1. 

Several observations can be made at this stage: 

(i) No spurious frequency is created if either the analytical inflow boundary conditions 
(2.23) - (2.24) qz the numerical one (2.25) are in characteristic form. 

If the analytical boundary conditions are in characteristic form, then ao = /3o = 0 and 
the system case is reduced to the scalar case. If the numerical inflow condition (2.26) 
is in characteristic form, then cr 0 = £o = 0 and the system again reduces to the scalar 
case. 


(3.4) 

(3.5) 
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(ii) On the other hand, if the inflow boundary conditions are not given in terms of char- 
acteristic variables, then there is a wealth of possible boundary treatments that lead 
each to different spurious solutions. In fact, for every frequency there exists a set of 
inflow boundary conditions that induces a spurious solution with that frequency. In 
particular, given an arbitrary zq = e*^; ^Owe get k 3 (z 0 ) from (2.30). We can then 
view (3.3) as an equation for a 0 cr 0 and Pq£q. Since these quantities are real, whereas 
the k s (z 0 ) are complex, Equation (3.3) is a system of two equations for a 0 a 0 and (3 0 e 0 , 
yielding a unique solution. The explicit form of the spurious solution (3.4) reveals the 
nature of such a solution. We summarize the results in the results in the following 
lemma: 


Lemma (3.1): The spurious solution (3-4) converges to zero for every fixed x = j Ax and 
t = nAi as Ax — > 0, A t —> 0. 


Proof: Since Ax — ■» 0 and x is fixed, then j — * oo. Note that |/c,j < 1, and therefore 
limj_* 0 k-J = 0. 


Notice that even though Lemma (3.1) indicates convergence of the total scheme, still on 
a fixed finite grid, the spurious solution 


j.n 

W 3 



may be confused with a real time periodic solution. In particular, if one of the /c 4 (e*^) is 
close to unity in magnitude, then the time periodic solution may show up in a large part 
of the spatial domain. Of course, if all the k’s are small in magnitude, the time periodic 
spurious solution will be confined to a very narrow boundary layer at the inflow boundary. 


To illustrate the above analysis, we chose to solve numerically equations (2.14) - (2.16) 
with Uoo = 1 and M ^ = .4. Thus, ai = —1.5; a 2 — 3.5; 03 = 1; = 2.5. For the initial 

conditions (2.18), we chose 

Ri(x, 0) = e~ x R 2 (x, 0) = e x R 3 (x, 0) = e 21 . (3.7) 

For the analytical boundary conditions (2.19) - (2.21), we chose 

a 0 = .582155; p 0 = -.6115; = 0; £i = 0 (3.8) 

and 

gi{t) = e _ajt - a 0 e° lt 

g 2 (t) = - Poe*' 1 (3.9) 

g 3 (t) = e ait . 

It is readily verified that under the above conditions the exact solution for the system is 

Ri{x,t) = e-( x " ait > 

R 2 (x,t) = e*-** 1 (3.10) 

R 3 {x,t) = e 2 ( I-ost ). 
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Thus the analytic solution decays to zero at steady state. The system was discretized by 
the Lax-Wendroff scheme (2.22) with the following parameters for the numerical boundary 
conditions (2.26) - (2.27) - (2.28) 

o=i = Pi = 0 
cr 0 = 1 e 0 = 2. 



Note that at the outflow boundary, x = 1, both the analytic and the numerical boundary 
conditions are in characteristic form, because in this section we are interested in the influence 
of inflow b.c.’s only. 

It can be shown that the inflow determinant condition (3.3) has for the above choice of 
parameters a solution of the form 

z Q = e 03 ‘. (3.12) 

The corresponding k,(z 0 ) obtained by the quadratic equation (2.30) satisfy 


|«!| = . 23035 |/c 2 | = .95334 |/c 3 | = .19015 (3.13) 

In the following graphs, we present w for several meshes 

1111 
X ~ 16’ 32’ 64’ 128 

At was taken to be O.lAx, corresponding to a CFL of 0.35 for our flow parameters. The 
results are sampled at x = .5, 0 < t < 9; thus 



jAx = .5 and 0 < nAt < 9. (3.15) 

At t = 9 the analytic solutions R, 9^ (s = 1, 2, 3) are of the order of 10 -6 . 

The numerical solution uy n versus time is presented for x = .5 and several meshes in 
figures: l(a,b,c,d) for = N = 16, 32, 64, 128. We do not present w^ ,n and w*' n because they 
converge even on coarse meshes ( N > 16) to the exact solution. This convergent behavior at 
x = .5 is due to the fact that |«i| and |/c 3 | are small and hence |/Ci|^ 2 and |k 3 |^ 2 are already 
of the same order (~ 10 -6 ) as the analytic solution. On the other hand, i u 2 ^J 2 displays a 
time periodic solution on each of the meshes presented. However, as expected, as the mesh 
is refined the spurious period is halved and the amplitude decreases. The solution converges 
with the decreasing mesh, but for every one of the fixed meshes used the solution at x = .5 
still exhibits the spurious frequency, because |/c 2 1^/ 2 is not yet small enough. This behavior 
at x = .5 is representative of the solution at other spatial locations. 


4. The Outflow Case 


In this section, we deal with the effects of non- characteristic formulation of the outflow 
boundary conditions in a similar manner to that of Section 3. 

It is easier to do the asymptotic ( N — ♦ oo) analysis in the outflow case by using new 
coefficients C’s. Defining 


C — 




(4.1) 


\C 3 J 
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we may rewrite (2.36) as follows: 


where 


A (ft) = 


A(/c)a + D{ji)c — 0 

(4.2) 

B(k)S+ E(ft)c= 0 

(4.3) 

D{li) = A(/i)A(/i) 

(4.4) 

Eiii) = fl(,.)A(,0 

(4.5) 

/ ft" +I 0 0 \ 


I 0 ft"!- 1 0 

(4.6) 


and the matrice A and B are defined in (2.35). 

Again, for N large, D(/i) and E(n) have small entries. In order to decouple the influence 
of the inflow conditions, we shall assume a= 0 and so the determinant condition becomes 

detl?(^) = 0 (4-7) 

or, more explicitly 

det E(n) = m (fj. 2 ~ 1)(^3 - 1) + ai^i/i 2 (Mi - 1)0*3 - 1) + - 1 )(t* 2 - 1) = 0. (4.8) 

A technical argument similar to the one leading to (3.4), (3.5) gives a solution of the form 

«»;•" = u ? 0 . . 

k(OI > i K 1 

The observation made there are also valid here, namely that a spurious solution (as defined 
there) may exist for non-characteristic outflow treatment; see observations (1) and (2) on 
pages 7-8. The numerical examples are again based on (2.14) - (2.16) with the same flow 
parameters, and with the same initial conditions given by (3.7). For the analytic boundary 
conditions, we took a 0 = fio = 0, o\ = l,£i = 2, and, from (2.23) - (2.25), 

<h = e-° at 
g 2 = e~ 2ast 

g 3 = e ait + (Tie -03 * + £ie -2ast . 

The analytic solution is again 

R 1 (x,t) = e-^~ ait \ i? 2 (M) = R 3 (x,t) = e 2 (*- a3t \ (4.9) 

The numerical boundary conditions (2.26) - (2.28) imposed on (2.22) were: 

a x = -4.668, fa = 2.09485, e 0 = a 0 = 0. (4.10) 

These values of ai and j3i were chosen so that we will have a case of a spurious frequency 
accompanied by at least one \l near to 1. The resulting fi’s and z 0 are: 


|/ii| = 1.00625, |;i a | = 2.0774, |/x 3 | = 1.3959, 2 0 = e al \ 


(4.11) 
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The numerical results for w *’ n are presented (again at x = .5) in Figures 2(a,b,c,d,e). 
and display a convergent behavior due to the smallness of \fJ’J N+3 \j=N/ 2 , { 3 = 2,3). On 
the other hand, is oscillatory on all the meshes we used, although the amplitude decays 
in accordance with the fJL\ N ^ 2 - law. 

We have thus demonstrated that spurious frequencies can be induced by improper (i.e., 
non-characteristic) treatment of boundary conditions, not only at the inflow boundary (see 
Ref. 1 for the physical case of a flow past a cylinder) but also at the outflow boundary. 

5. The Finite Domain Case 

In this section, we report on numerical experiments in which both boundaries (x = 0 and 
x = 1) are treated in non-characteristic fashion. We chose the same boundary conditions 
reported on in the previous two sections, but now they are applied simultaneously. Recall 
that when only the inflow was non-characteristic, the solution exhibited spurious oscillations 
but the solution converged as Ax — » 0, At — » 0. The outflow case on its own exhibited the 
same behavior. We mention also that the G-K-S theory [2], for the stable case (see Lemma 

(2.2) ), predicts that if the semi-infinite cases (i.e., pure inflow or pure outflow) are each 
stable separately, then also the finite domain case will be stable and convergent with mesh, 
though perturbation might grow in time. In the numerical experimentation reported herein, 
it will be seen that in our special case (|zo| = 1 for each of the semi-infinite problems) the 
solutions seem to grow exponentially in time and to converge, but non-uniformly, with mesh 
size. 

The numerical approximation (2.22) (for s = 1,2,3) was used to solve (2.14) - (2.16) 
with the boundary conditions now given by: 

a 0 = 0.582155 (3 0 = -0.6115 a 0 = l e 0 = 2 (5.1) 

ai = -4.668 A = 2.09485 a x = 1 e x = 2 (5.2) 

The boundary parameters (4.1) were taken from the “pure” inflow problem and those in 

(4.2) from the “pure” outflow problem. 

We show the temporal behavior of tu*‘ n at x = .5 in Figures 3(a,b,c,d,e,f). Examination 
of these figures shows that unlike the case of “pure” inflow and outflow, where there were 
oscillations but they did not grow in time in the present “finite domain” boundary treatment 
the oscillations grow temporally. w 2 ’ n and w 2,n display the same behavior and are not shown 
here. Another difference between the “pure” case and the present finite-domain case has 
to do with behavior as Ax is being decreased. In the previous sections, we saw that the 
amplitude of the oscillation decreased uniformly with the mesh size. In the present case, 
there is an increase in amplitude as Ax is decreased from N = 16 to N = 128. Afterwards, 
for N = 256 and 512 we see a decreased amplitude; and although the perturbing oscillations 
are still large compare to 10 -6 (and growing temporally) there does seem to be a beginning 
of convergence. 

We should remark here that this type of behavior of the numerical approximation to the 
solution can be very unsettling for practitioners using practical codes, since multidimensional 
mesh sizes seem likely not to be in the asymptotic range in the near future. 
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Figure la. vs. t < 9 (Inflow, N = 16) 
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Figure lb. w]^ 3 vs. t < 9 (Inflow, N = 32) 





Figure Id. w ^" 2 vs. t < 9 (Inflow, N = 128) 
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Figure 2a. u>jy " 3 vs. t < 9 (Outflow, N = 16) 
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